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arising from the interaction between a pair of vorticity-propagating waves. To illustrate and elucidate the physical ef¬ 
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instability. We find that the non-Boussinesq effects introduces a preference for the direction of wave propagation that 
depends on the sign of the shear in the region where waves are supported. This in turn affects phase-locking of waves 
that is crucial for the mechanistic interpretation for shear instability, and is seen here to have an inherent tendency for 
stabilisation. 
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I. INTRODUCTION 


An approximation that one often makes when studying the dynamics of stratified fluids is the Boussinesq approximation^*. 
One assumes that the variation of density about a background reference density p' / po is small, and thus we may neglect inertial 
effects associated with such terms except when it is multiplied by the gravitational acceleration g, i.e., buoyancy effects dominate. 
This assumption of small density deviation is well satisfied in the ocean and remains useful for studying certain atmospheric 
flows. The Boussinesq equation has and still remains a useful model for investigating a variety of fluid dynamical phenomena in 
geophysical systems, such as convection^, wave-mean flow interaction^ and shear instabilities^, the last of which will be our 
principal focus here. 

To study shear instabilities, one often makes a further simplifying assumption by employing ‘defects’ in the velocity and/or 
density profile (i.e., piecewise-constant/linear profiles) as a model for sharp gradients in the basic state. Such an assumption is 
useful for the study of the onset of instabilities for several reasons: the resulting dispersion relation often reduces to a low order 
algebraic equation, for which analytical as well as asymptotic solutions exist; such solutions are often the leading asymptotic 
solution for general smooth profiles in the long-wave limit^; there is a mechanistic interpretation for the instability, seen as the 
constructive interference of vorticity propagating waves travelling counter to the background flow^. The use of defects has life 
beyond linear theory, allowing the derivation of reduced models via matched asymptotic methods to investigate the nonlinear 
development and saturation of shear instabilities!^. Since the use of defects is as a model for sharp gradients in the basic state, 
one can ask whether it might be more appropriate to study flow instabilities in the presence of sharp density gradients without 
the Boussinesq approximation, since the assumption of small density variation may no longer hold. To this end, there have 
been several works studying shear instabilities beyond the Boussinesq approximation over the years, using smooth profiles but 
with a density that has a small scale height* 1 , and classic profiles with defects irt^^i. However, these aforementioned works 
in the non-Boussinesq setting focuses on solving the modified Taylor-Goldstein equation to investigate the property of growth 
rates with increasing deviation from the Boussinesq regime (which will be seen to be measured by a Froude number), without 
necessarily providing a physical reason of what causes the modifications to the instability characteristics. Our work here aims 
to complement these previous works by investigating the mechanistic modifications to the underlying wave dynamics by non- 
Boussinesq effects, and how this affects the mechanistic interpretation of the instability accordingly. We provide mathematical 
details and physical schematics on how the part of the baroclinic torque neglected by the Boussinesq approximation generates 
vorticity anomalies; how this affects wave propagation and interaction is illustrated for increasingly more complex examples. An 
instability problem where the cause of instability is strongly affected by the non-Boussinesq term is then presented and analysed 
accordingly. 

The layout of the document is as follows. In Section|I]] we formulate the problem in terms of the vorticity, displacement and 
pressure, to relate the generation of vorticity anomalies by the Boussinesq and non-Boussinesq effects. The dynamics of waves 
supported on a uniform background is investigated and rationalised in Section!]]]! to illustrate some of the possible effects due to 
the non-Boussinesq term. In Section [IV] we consider a more complex example where waves are now supported on defects, and 
rationalise also the changes induced by the non-Boussinesq term. In Section[V] a slightly simpler version of the Taylor-Caulfield 
instability-^^ in the non-Boussinesq regime is investigated and analysed accordingly. This ties together the modification to 
the wave dynamics and the instability properties resulting from the action-at-a-distance interaction between non-Buossinesq 
interfacial waves. We conclude and discuss our results in SectionlVIl 


II. MATHEMATICAL FORMULATION 


We assume a two-dimensional, inviscid incompressible flow in the (x. z) plane, with governing equations 


Du 

~Dt 


= -ge z 




(1) 


where u = (u, 0, w) is the velocity, g is the gravitational acceleration, p is the kinematic pressure, p is the density, and D/Dt = 
d/dt + u ■ V is the material derivative. The last equation for density comes from assuming V ■ u - 0. Defining q = 
dw/dx — du/dz (note that the vorticity component in the y direction is uj y = —q), the q equation is given by 


Dq 

Dt, 
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1 / dp dp dp dp\ 
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How both terms of the baroclinic torque generate vorticity anomalies is illustrated in the schematic depicted in Figure [I] with 
the details in the caption. 

Suppose we take a basic state U(z), with Q(z) = —dU jdz. Then it follows from the ^-component of the momentum equation 
that we may take p = p(z). This leads to the basic state satisfying the hydrostatic balance 


dp 

dz 


= -pg, 


(3) 
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FIG. 1. Vorticity generation via baroclinic torque associated with p and p anomalies. We show here a case where ( dp/dz){dp/dx ) > 0, a 
case resembling hydrostatic balance, and a case where — (dp/dx)(dp/dz) > 0. Pressure gradients results in a tendency in flow acceleration 
flow. However, the acceleration is inversely proportional to the density of the fluid parcel, thus it results in shear and hence vorticity anomalies. 


and we consider a basic state p = p(z). With this, we observe that the terms in the square brackets above will be O(e) once we 
linearise, and so the contributions at O(e) from p~ 2 is (p) -2 . A linearisation of the vorticity equation then results in 


. tt ^ i dQ 

— + U — ) q = — w —— 


dt 


dx 


dz 


1 /_ dp dp dp\ 
p 2 dx + dz dx) 


(4) 


where the quantities with no overbars are perturbation quantities. Linearising the continuity equation results in the system of 
equations given by 


/ = _ dQ g dp 1 dp dp 

\dt dx) q W dz pdx -fdzdx ’ 





(5) 


upon using the divergence-free condition on the perturbation velocity field. Note that the last term on the right hand side of the 
linearised q equation is absent in the Boussinesq limit. 

With the vertical perturbation displacement defined as 


d_ 

dt 


U 


-V- 


dx J 


( 6 ) 


an integration yields the identity 



(7) 


where only advective effects are considered. Since the velocity field is assumed to be non-divergent, we may define a stream- 
function such that 


u = — 


dp 

a? 



( 8 ) 


and this results in the identity V 2 ip = q. Defining the Buoyancy frequency to be N 2 = —(g/p)(dp/dz), the system of equations 
© becomes 
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(9) 


The (d/dx){p/p) term is the correction that is absent in the Boussinesq regime. Here, tb may be formally inverted from q via a 
Green’s function, which depends on the chosen domain and boundary conditions, so in theory we have a formulation in terms of 
q and tb, once we substitute for p in some way. The pressure p will be seen to be related to tp and thus q via a substitution from 
the momentum equation. 

We observe that there are three dynamical regimes: 
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1. A barotropic regime where ( = p/ ( gp ). For barotropic flow, p = p{p), so that 


Dp{p) = Ap_Dp , 1Q , 

Dt d p Dt 

With incompressibility, Dp/Dt = 0, and so p = —(/{dp/ dz). Upon using the hydrostatic balance relation ©, we recover 
the identity £ = p/{gp). In this regime, the linearised baroclinic torque cancels out exactly, consistent with the assumption 
of barotropicity. 

2. The Boussinesq regime where ( p/{g7>)- In this regime the fluid parcel adjusts its pressure distribution to the sur¬ 

rounding environment on a fast enough time-scale such that pressure effects may be neglected, and buoyancy effects 
dominate. 

3. The case where / <C p/{g~p). In this case, the non-Boussinesq effects outweigh the buoyancy effects and become the 
dominant player in the dynamics. 


We will now consider related examples of increasing complexity to see how the extra non-Boussinesq term in the linearised 
baroclinic torque influences the dynamics. 


III. BASIC WAVE DYNAMICS 


We consider first the case where waves are supported on a uniform background flow with N = constant. Without loss of 
generality, we take U = 0, and the governing equations © reduces to 


d/l = ( c _ }_p\ dC_ = dtp_ 

dt dx \ gp / dt dx 


( 11 ) 


To substitute forp, we turn to the momentum equation du/dt = —{l/p){dp/dx). Substituting for p in equation (fTTb . taking 
another time-derivative of the vorticity equation and substituting for d(/dt results in 


dt 2 



N 2 dip\ 
g dz) 



( 12 ) 


This is the Taylor-Goldstein equation for this simplified case. 

Now, taking p = poe ~ z / H , dp/dz = —~p/H, so N 2 = —{g/p) {dp/dz) = g/H, where H is a density scale height. Equation 
£□ becomes 


With modal solutions of the form if) 



t/){z)e ife (*-ct), 



-i®. = 


H dz 


dx 2 


we obtain the dispersion relation 
1 d 


H dz 


-k z )'ip = —Ni/j. 


(13) 


(14) 


We see that solutions of the form /p ~ e imz e z/( 2 H) sa ti s fl es \u \ 2 ^ e z / H and p ~ e Z ^ H , so p\u \ 2 < oo. Substituting this form 
of solution into equation (IT4l) . the imaginary parts cancel out exactly, and the resulting dispersion relation is given by 


c 


2 


N 2 

m 2 + k 2 + l/(4fT 2 ) 


(15) 


In the Boussinesq limit, H —> oo, and we recover the usual dispersion relation for gravity waves in a non-rotating system^. For 
H < oo, the phase speed of the waves reduced via the l/(4fT 2 ) term. This is akin to how the scale height affects acoustic- 
gravity waves^, and similar to the way in which the existence of a finite Rossby deformation radius attenuates the phase speed 
of Rossby wavesaS. 

A physical reason for this reduced phase speed may be rationalised via the changes to vorticity anomalies generated by the 
corresponding baroclinic torques. For simplicity, we consider the case with m = 0, and, for completeness, we consider the 
Boussinesq limit first. The vorticity generation comes from the {dp/dz){dp/dx) = —g~p{dp/dx) term; since g~p > 0, the sign 
of the resulting vorticity anomalies is correlated with the sign of ~{dp/dx). The direction of the wave propagation is dependent 
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FIG. 2. Schematic for gravity waves propagation in the Boussinesq limit, (a) Direction of propagation depends on how ( and q are correlated 
at the peaks and troughs. Shown here is the case where ( ~ q, and results in c > 0; with £ ~ —q, the opposite is true. (6) Vorticity anomalies 
at the nodes resulting from the baroclinic torque associated with the Boussinesq term, which depends on the £ configuration and not on the 
direction of wave propagation. 


on how (is correlated with q; £ ~ q gives rightward propagating waves with c > 0, as in Figure[2] «), and vice versaai. However, 
regardless of direction of wave propagation, equation £7} indicates that, for stable stratification, the £ distribution sketched in 
Figure |26) results in dp/dx < 0, and thus this results in positive vorticity anomaly at the node of the wave. Another way of 
thinking about it is that the peaks of the wave has the tendency to descend whilst the troughs wants to rise, therefore the resulting 
movement in this case is anti-clockwise, and is thus a positive vorticity anomaly as in FigureQJn). 

In the more general case with the non-Boussinesq term we also need to work out the distribution of p and see how this modifies 
the scenario depicted in Figure|2 From the x-momentum equation, for a right-going wave (c > 0), observing that ip ~ e z /(' 2lr >_ 
we have 


dtp p 

— cu = c— = - = 
dz p 



-ip ~ q, 


(16) 


and thus p ~ q&. This scenario is depicted in Figure|2a). The resulting p distribution leads to dp/dx < 0 at the nodes, and thus 
the vorticity anomaly generated is related to —{dp/dz){dp/dx) < 0, resulting in a negative vorticity anomaly at the node as in 
FigureQ] b). So the correction torque results in vorticity anomalies that is of the opposite sign to the one generated by Boussinesq 
term depicted in Figure[26). This may be seen to reduce the wave propagation speed since the speed is related to the magnitude 
of the vorticity anomaly generation at the nodes^I. The magnitude of the pressure anomalies and thus the resulting vorticity 
anomalies are related by the size of H , which in this instance measures the degree of deviation away from the Boussinesq limit. 
The same line of thought may be applied to the c < 0 case, which results in this q ~ —p, and the corresponding scenario is 
illustrated in Figure[26). Again, the sign of the resulting vorticity anomaly is seen to be opposite to the one given in Figure[2&). 
Thus, in this setting, the physical picture is that the baroclinic torque associated with the non-Boussinesq effects reduces the 
wave propagation speed via generation of opposite signed vorticity anomalies to the ones generated by the Boussinesq term at 
the nodes. This reduction is symmetric in magnitude for waves propagating in either direction, since the pressure distribution 
associated with the waves with c > 0 and c < 0 depends only on the distribution of </ (as p ~ c(, which may be seen when 

combining equations fill and IT6l>. We will see in the next section how the background shear affects the propagation of right and 
left going waves in an asymmetric way. 


IY. EDGE WAVE DYNAMICS 

Suppose now our background profiles are piecewise-continuous, so that dQ/dz and dp/dz (and so N 2 {z )) are defects of the 
form 


dQ_ 

dz 


A QS(z — h), 


N 2 {z) = A N 2 S{z-h). 


(17) 
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FIG. 3. Schematic of the torque generated by the non-Boussinesq contribution for gravity waves with m = 0. The imposed wave structure is 
as in Figure [2] the left and right panels depict a right and leftward propagation wave respectively, and the shading is associated with contours 
of ip = const. As a result of the vorticity distribution at the peaks and troughs, we have the corresponding pressure distributions which in turn 
generates vorticity anomalies at the wave nodes; these are seen to be the opposite sign as the one associated with the Boussnesq term shown in 
Figure |2} 6). 


Then it may be seen that solutions of the form 


q = qe ik ( x ~ ct) 6{z - h) (18) 

are consistent solutions of (0 since there is no vorticity generation away from the location of the defect at z = h. Taking also 
modal solutions of ip, ( and p, we have 

AN 2 n 

(U - c)q = -AQip - AN 2 ( + -£ [U — c)£ = ip, (19) 

a p 

where all the relevant terms are to be evaluated at z = h. In a domain that is unbounded in z, ip is related to q via a Green’s 
function 


ip(h) 


! 

2k 


q(z')e-^ h - z ' l dz' 


a_ 

2 fc’ 


( 20 ) 


and it remains to relate p to the prognostic variables q and 

We now wish to substitute p/~p for u = —dip/dz by making use of the ^-momentum equation. Generically, ip is not differen¬ 
tiable at z = h, however, u changes sign when z = h is crossed, hence, physically, for a wave supported on z = h, there can be 
no self-induced u for a wave-like solution, and thus u = 0 at 2 = h. With this, the x-momentum equation becomes in this case 


o=-f?+0^ 

ox p ox 



( 21 ) 


This means that, in the absence of shear, the pressure perturbation of the interfacial wave is zero, and the Boussinesq approxi¬ 
mation holds exactly in the linearised baroclinc torque. 

Substituting for ip and p/~p in (fl9t . we obtain 


(U -c)q=[AQ- ± - AN 2 C, 




( 22 ) 


The eigenstructure and dispersion relations are thus given by 


= 2fc(c ± - [/X , ( c± — U) = — 


AQ-(Q/g)AN 2 \ 
4 k ) 


i 


( AQ — (Q/g)AN 2 
v 4/c 


AN 2 
2 k 


(23) 


6 















p-Ap 



FIG. 4. Schematic of the torque generated by the non-Boussinesq contribution for edge waves. With positive shear, we have Q < 0, and thus 
—ip q r*j p from GD The choice of shear fixes the relation between q and p at the peaks and troughs. With this, the associated pressure 
distribution generates vorticity anomalies at the nodes, introducing an asymmetry to the direction of wave propagation, in this case to leftward 
propagation. 


The coefficient ( Q/g)AN 2 measures the deviation away from the Boussinesq limit. When AN 2 = 0, we recover Rossby waves, 
while for Q = AQ = 0, we recover the gravity waves in the Boussinesq regime?-. The plus and minus branch are the branches 
where the appropriate sign is taken. 

We make the observation that, for ( Q/g)AN 2 (AN 2 , AQ), one of the branches vanish and so the edge waves become 

uni-directional, with the preference direction dependent on the sign of Q. This is like the case for the propagation of Rossby 
waves, where in the absence of stratification, is uni-directional and depends on the sign of AQ. This preference for the direction 
of propagation is caused by the extra contribution to the baroclinic torque from the non-Boussinesq term. To see how this 
operates, we consider the case where we have a positive shear with Q < 0. With this chioce, q ~ —ip ~ p from (12Tb and this 
fixes the p and q relation at the peaks and troughs, as displayed in Figure [4] Then we may consider both the case where q ~ £ 
(for rightward propagating waves) and q ~ — £ (for leftward propagating waves). Taking into account Figure[2j6), it may be seen 
that the resulting pressure anomalies results in vorticity anomalies at the nodes that is of the opposite sign to the base case for the 
rightward propagating wave (base case of Figure 0 a)), i.e., counteracts propagation; on the other hand, the vorticity anomalies 
at the nodes for the leftward propagation wave is the same sign as the base case (of Figure 06)), i.e., reinforcing propagation. 
Thus there is a leftward preference for wave propagation due to the non-Boussinesq contribution with positive shear; this is seen 
to be consistent with the dispersion relation given in (|23] >. 

In this scenario with edge waves, the p and q relation is fixed by the sign of the shear, where as in the Section [II] for neutral 
waves supported in a flow with no shear, the distribution depends on the direction of wave propagation. The correction to the 
baroclinic torque acts to counteract wave propagation in a symmetric way for the neutral wave case, whilst the presence of a 
shear introduces a preferred direction for wave propagation. 


V. NON-BOUSSINESQ TAYLOR-CAULFIELD INSTABILITY 

One mechanistic interpretation for the onset of shear instabilities is via the constructive interference of counter-propagating 
waves. Waves that propagate vorticity anomalies may become phase-locked with each other via the advection by the background 
flow and action-at-a-distance of the nonlocal velocity field induced by local vorticity anomalies. With phase-locking, depending 
on the phase shifts, these waves may amplify each other and lead to instabilit y 7,18 ! 21 ' 23 . 

Since one of the key ingredients for this interpretation is counter-propagation, our hypothesis with the previous section in 
mind is that, when the correction term (as measured by with relation to {Q/g)AN 2 ) becomes significant, instabilities reduce in 
growth rates and eventually switch off because the waves can no longer phase-lock as they become increasingly uni-directional. 
This suggests a physical interpretation to the work of Barros & ChoUi, who find that a large shear across the interfaces plays a 
stabilising role, which is perhaps somewhat counter-intuitive as the shear is normally seen as a source of instability. To test this 
hypothesis, we consider a simplified form of the the Taylor-Caulfield proble m 10,15 ' 16 ' 18 ! 24 , where the basic state is essentially 
given by 


U(z) = A z, N 2 (z) = AN 2 [S(z -h) + S(z + h)}, 


(24) 
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with the 5-functions in N 2 (z) coming from the choice that p = po + |< 5 p|[l — H(z — h) — H(z + /i)], H(z) the Heaviside function, 
h are the locations of the defects, and the imposed density is a staircase-like profile. With this, the instability comes from the 
interaction of two interfacial gravity waves located on the defect a^& — 

We proceed to non-dimensionalise the equations. By scaling with respect to To = A -1 and Lq = h, and taking modal 
solutions as in (fl8l . it may be seen that the dimensional equations (0 becomes (noting that (5-functions have dimensions L^ 1 
and that dQ/dz = 0 here) 


(±1 — c)<h,2 = —R ^Ci ,2 — F 2 (±1 - c)Ct ,2 = Vh, 2 , (25) 

where the equations are evaluated at 2 = ±1 for subscript 1 and 2 respectively, and all quantities are non-dimensinonal. The 
non-dimensional parameters in this case are 


~ AN 2 _ h 2 A 2 


(26) 


The Richardson number R measures the strength of the stratification. The Froude number F is given by the square of the mean 
shear velocity scaled by the Boussinesq gravity wave speed. Since the presence of shear allows the non-Boussniesq baroclinic 
term to operate, it measures the deviation from the Boussinesq limit. In the limit F —» 0, we recover the Boussinesq limit where 
solutions to the problem as stated are knowni&. With this rescaling, the edge wave structure (l23l > associated with this set up is is 
given by 


qf 2 = 2k(c ± - U) ia Ci t2 i 



R 

2k' 



(27) 


It remains to relate ipi i2 and (p/p) i i2 to qi^ and £ 1 , 2 - 
unbounded domain as in < [20b . except here we have^ ! 


First, ip 1,2 may be related to q\ 2 via the Green’s function in an 


Vh,2 — ~2fc(7i,2 + 92,ie 2fe )- (28) 

Note that we have a term with a flipped subscript to denote the interaction induced by anomalies on the other interface, with the 
exponential factor representing the decay of interaction strength. For (p/p) 1 , 2 , we again make use of the ^-momentum equation, 
which is, in this setting and with u = —di/t/dz. 


~ (±1 - c) 


dip 1,2 
dz 



(29) 


The physical argument here is that there should be no self-induced u on an interface but there may be an induced u from the 
other interface. Since a positive vorticity anomaly induces a positive u below and negative u above it (and vice-versa for negative 
vorticity anomalies), we obtain 


dip _ f-i / q(z')e fe(z z,) dz', z > z', dipi^ _ 1 „ _ 2 k 

dz [ + ! / q(z')e +k ( z ~ z ') dz', z < z 1 , ^ dz 

Substituting the above into (f25l> . we obtain the governing system of equations 


(30) 


(il — c)(gi ! 2 T 7^?2,ie 2k ) — — 7?Cl,2 + 7 ( 91,2 + 92,te 2k ), (±1 — c)Ci,2 — — ^( 91,2 + <72,ie 2k )- (31) 

In matrix form, this is 


/ 

1 — c — 7 

R 

—'ye 2fe [l + (1 — c)k] 

° ^ 

/< 7 t\ 


l/(2*0 

(1-c) 

e ~ 2k /(2 k) 

0 

Cl 


-7e _2fe [l + (1 + c)k\ 

0 

-(1 + C + 7) 

R 

V 2 J 

V 

e ~ 2k /(2 k) 

0 

l/(2fc) 

-(1 + c)/ 

\C2/ 
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( 32 ) 













FIG. 5. Contour of growth rate kd in (k, R) space at various values of F; all the contours are at fixed levels at 0.01 spacing for all four panels. 
The blue dashed contours are where the resonance condition ( 134b is satisfied. 


and this yields the dispersion relation 


(1 + 7 2 e~ 4h k 2 )c 4 + 2 7 c J - 2 + (R/k) - 7 2 (1 - e“ 4fc (l + 2 k + 2k 2 )) 


2 7 (1 - f?e _4fe ) + 7 (i?//c)(l - e _4fe ) 


1 “ I + S (1 “ + ^( e ' 4fe ( 1 + fc ) 2 - !) 


= 0 . 


(33) 


When F = 0, 7 = 0, and (l32l) as well as (f33l) reduce to previously known forms in the Boussinesq limifc^. 

The dispersion relation ( l33t may be solved numerically to obtain the four roots and this was done using the MATLAB 
command roots. We first show in Figure0the contours of the growth rates over (k. R) space at several values of F. Figure[5ja) 
is exactly the solution in the Boussinesq limit for which an analytic expression for the solution is available^. As we increase F, 
the growth rates reduces in Figure[5j &), notably around the region of maximum growth. As F is increased further, the maximum 
growth rate decreases, and the region of instability shrinks towards the small k region, as seen in Figure [^c, d). 

Sometimes it is useful to show the locations where the resonance condition is satisfied 7 . These are the locations where the 
counter-propagating edge waves have matching phase speeds, taking into account advection by the background flow. From 
equation ( fTTI i. these are the values of k where 


c i 


(R 


> F, k) — 1 + 




and 


c+{R,F,k) = -l + 




(34) 


are equal. If interacting counter-propagating edge waves contribute the most to the dynamics, then the location where the 
resonance condition is satified should be near to the location of optimal growth; otherwise, it shows that other dynamics (e.g. 
pro-propagating modes, critical layers) are important. It also gives an indication of where in parameter space the interaction 
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FIG. 6. Line graphs of the solution branches of the dispersion relation 1331 . with blue solid lines denoting c r , and red dashed line denoting 
20 ka (factor of 20 to emphasise the instability region), at different values of R and F. The rows are at different values of R whilst the columns 
are at different values of F (see diagram). 


required for instability may be expected. Locations of these are shown as dashed contours in Figure [5] and we see these show 
reasonable correlation to the locations of largest growth. However, we notice that, in the larger F cases, even though we may 
have edge waves with matching phase speeds, this does not necessary indicate instability, since the resonance condition does not 
take into account the mutual interaction. 

It is perhaps informative to see the behaviour of the individual solution branches. In Figure [6] we show the phase speeds c r 
(solid blue) and the (magnified) growth rate 20 kci (dashed red) for several values of R and F. The first column shown in panels 
(a, e, i, m) is the Boussinesq case where there is no preference for the direction of wave propagation for this choice of basic 
state. Focusing on the phase speed, the outer two branches are the neutral waves, while the inner branches represent the stable 
and unstable branch, occurring in conjugate pairs. As F increases, the degree of asymmetry increases, with a preference for 
leftward propagating waves, which is consistent with the result in Section[IV] We also note that although the branches appear to 
cross, there is only instability when a pitchfork-like (rather than a transcritical one as in panels ( d , h, l )) bifurcation occurs. 

To further quantify the asymmetry between the leftward and rightward propagating waves, we wish to obtain the instability in 
terms of the left and rightward propagating modes Unlike the previous formulations where transformation matrices were 
present^ 1 :, the complication here is from the (p/p) 1,2 term in equation (l25l >. This contributes a ( d/dt + Ud/dx)u term, which 
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means we can no longer write the problem in the form dC/dt = A£ in a simple way, and the transformation matrix acting on 
A becomes complicated. In principle, since everything is linear, an alternative approach that one could take is to work out how 
the individual terms in the governing equation 031) should look like, and the equations for Cy/. 2 should have on the right hand 
side the interaction terms written in terms of the appropriate contributions from the terms in (ESI . We may postulate for example 
that the equation for say should only be affected by all variables not including i.e., the governing equations without the 
modal solution assumption should be of the form 


(±1 - c)qf 2 = -Ri ,2 



(±1 - c)Ct,2 = 


dip* 

~dx’ 


(35) 


where 


r = i>t,2 +{i’ll + ^2,1) e ~ 2k > 


(36) 


and (p/p)* is to be defined analogously. With this, we may substitute accordingly noting that: (i) q/p 2 — > 2fc(c ± — I/) 1 , 2 Cm: y i a 
the eigenstructure (ETI : (ii) we take ipf 2 = —qf 2 /2k, and again may be written in terms of (p/p 2 via the eigenstructure; (iii) £±2 
is a local variable so we leave it as is; (iv) some care needs to be taken for the (p/~p)* term, but we essentially use the definition 
that 


P(z) = 1 _ 10, z = z b 

P \(c(z b )-U(z))(dip/dz), z^z b , 


and these may be then be written in terms of £f 2 via appropriate substitutions. 

The resulting manipulations are quite unwieldy due to the large number of terms and we shall not present them here. Instead, 
we may achieve the same goal by decomposing the resulting unstable modes into its normal modes. Since we already have c 
from the calculations, one way to do this is to write (1321 as 
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and solve for c and (</\. £ 1 , c/ 2 , £ 2 ) using the MATLAB eig command, but returning a solution only |c — cl < 10 15 (this 
condition is satisfied for all unstable solutions computed here). To then transform the resulting solution into £f. 2 , we make use 
of the wave structure (1271 . so that the unstable mode in terms of normal modes is given by 21 




cr 


c 2 + 


\C 2 / 

\ 


2 k(c~ - U) 1 
1 
0 
0 


0 0 \ 1 

0 _ 0 

2 k(c+ - U) 2 2ft(cr - U) 2 

1 1 / 



(39) 


The plus and minus superscripts denote the rightward and leftward propagating modes, and it is primarily the interaction between 
the two counter-propagating modes and ( 2 that leads to instability, with the pro-propagating modes ££ and ( 2 that modify 
the interactions accordingly. 

With 




± 

1,2 


A; 


±„ie 


± 

1,2 


(40) 


we define Af > 0 and ef 2 € (—tt, tt] to be the (real) amplitude and phase of the respective modes. With this, we show in 
Figure [7] the (normalised) phase difference Ae/n between the counter-propagating modes and the ratio of the amplitude of the 
rightward-propagating waves and the leftward-propagating waves r, respectively given by 


Ae = e% -e 1 , 


A+ + A+ 

T = — -—. 

Ai + A 2 


(41) 


Starting first with the phase difference, since we defined it using the displacement rather than vorticity—, it may be seen that 
—tt < Ae < 0 is the unstable regime, and with — 7 r /2 < Ae < 0, we are in the ‘hindering’ regime where the configuration 
is such that the counter-propagating waves hinder each other’s propagation to achieve phase-locking, typical of fast waves (see 
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FIG. 7. Diagnostics from the normal modes. (a,c,e,g) shows the normalised phase difference Ae/ir, where —0.5 < Ae/ir < 0 is the 
‘unstable hindering’ regime, and —1 < Ae/ir < —0.5 is the ‘unstable helping’ regime. ( b , d, f, h ) shows the ratio of the total amplitude of 
rightward-propagating modes and the total amplitude of leftward-propagating modes. 


also Figure [8]here). We see in panels (a, d, g, j) this occurs for waves at higher values of R and small k, which is consistent with 
the dispersion relation d27} . where faster waves occurs for larger R and smaller k. The reverse is true when —n < Ae < —ir/% 
and we are in the ‘helping’ regime. This explains why the locations where the resonance condition (134} is satisfied does not 
necessarily correspond to the location of largest growth rate. As indicated from the Green’s function (l20l > and < 128 b . the interaction 
strength between the waves increases as k decreases, but then so does the counter-propagation speed from ( 123} and {27}. Hence, 
the gravest mode is obtained in growing, hindering configurations, which is a generic result that applies to barotropic and 
baroclinic instabilitie s 26 ' 27 . 

For the ratio r as defined in ( ITO . we make the observation that, for the Boussinesq limit where F = 0, there is no preference 
for direction of wave-propagation, so the value of r should be equal to 1 over the parameter space, which is what we see in 
panel ( b ). For non-zero F, there is a preference for leftward-propagation, so the value of r is less than 1 and decreases in size 
as F increases, which is what we observe in panels (d, f, h). We make the observation that the asymmetry is less strong for 
long-waves, indicating the non-Boussinesq effect appears to have a stronger effect on short-waves. This is perhaps consistent 
with the expectation that we expect buoyancy effects to remain dominant for large-scale motions, and non-Boussinesq term 
affect small-scale motion more substantially. Notice that the r does not need to vanish for the instability to switch off; waves 
being unidirectional is a sufficient but not necessary condition for phase-locking, and the ability to phase-lock may disappear 
before waves become unidirectional. 

In the work of iRabinovich et a/3 ^ in the Boussinesq regime, it was argued that, for phase-locking, the pro-propagating mode 
on one flank should be in anti-phase and smaller by a factor of x with the counter-propagating mode on the other flank. The 
picture is likely to be somewhat more complicated here in the non-Boussinesq regime. There is now a preference for the 
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FIG. 8. A regime diagram and for phase differences and the likely physical scenario of the resulting instability in the non-Boussinesq regime. 
The regime diagram is defined using displacement with phase difference as defined in equation G3. The wave configuration between the two 
counter-propagating modes is in a growing hindering regime here, with the top wave travelling to the left and bottom wave travelling to the 
right. The top wave is of large amplitude, as indicated by the larger labels and arrows, consistent with the diagnostics as in Figure[7] Although 
the waves are phase-locked, the leftward propagating wave dominates and imparts a leftward propagation to the resulting instability, consistent 
with the observations in Figure [6] 


direction of travel, and thus may not be (and is generically not) equal to X '2 — A 7 /A^ • These diagnostics do 

not tell us anything overly meaningful, and a presentation of the associated results has been omitted. For completeness, the 
ratio between the two counter-propagating modes xo = /^i shows that they are mostly comparable in magnitude over 

the unstable region, with a slight preference towards the leftward counter-propagating mode A~. As we have seen before, 
the observed instabilities are no longer stationary modes, and thus a clear picture as in iRabinovich et al . I is unlikely to hold 
in this case. A likely physical scenario for instability is that phase-locking is still achieved, but, as we may expect from the 
hypothesis, since there is a preference for leftward propagation, the leftward travelling counter-propagating mode dominates 
and imparts a leftward propagation to the resulting instability, which is consistent with the negative values of c r observed in 
Figure 0 A schematic of the resulting interaction between counter-propagating waves is shown in Figure [8] and we expect this 
to be the fundamental component in driving the instability, with the pro-propagating modes modifying the interaction in a more 
complicated manner. 


VI. CONCLUSION AND DISCUSSION 

In this article, we investigated how the portion of the baroclinic torque that is neglected by the Boussinesq approximation 
affects wave propagation, and how this in turn affects the mechanistic interpretation for shear instability. Increasingly complex 
examples were considered and rationalised, and it was found that the dynamics depends on the relation between pressure and 
vorticity anomalies. In Section [III] we observed that, for neutral waves supported on a uniform background flow, the non- 
Boussinesq term acts to reduce the wave speed via generating vorticity anomalies that counteract what would otherwise be 
generated by the Boussinesq term, in a symmetric way that depends only on the direction of wave propagation. In Section ITVl 
the introduction of a background shear fixes this degree of freedom between the pressure and vorticity anomalies and introduces 
an asymmetry for direction of wave propagation, which is to the left for positive shear (Q < 0). In Section [V] a simplified 
version of the Taylor-Caulfield problem was investigated and analysed. With positive shear, the hypothesis was that, since 
there is a preference for leftward propagation of waves, as we increase the non-Boussinesq effect as measured by the Froude 
number F, the waves should become increasingly uni-directional. With this, phase-locking becomes harder to achieve, and thus 
increasing F reduces the region of instability and the maximum growth rates. This was indeed found to be the case via plots of 
the maximum growth rate in Figure [5] and the values of the ratio of the total rightward propagating waves to the total leftward 
propagating waves r shown in Figure[7] the latter obtained by a decomposition of the unstable modes into its left and rightward 
propagating constituents via the dispersion relation (l27ll . 

These results are in general agreement with the previous works on shear instability in non-Boussinesq systems^ - — even if 
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their precise set up is not identical to ours. In particular, the observation that increasing the shear (i.e., the value of F) stabilises 
the instability^ is in agreement of our results here, with the reason being that the magnitude of the shear increases the degree 
of asymmetry for wave propagation, which in turn affects phase-locking properties. Furthermore, the work of Barros & Choi^ 
analysed the non-Boussinesq effect when it is combined with the effect of confinement by boundaries. As shown in some 
previous work s 28 ' 29 , the reduction in growth owing to confinement can also be explained in terms of wave interaction, since 
mirror image waves that are in anti-phase with the counter-propagating waves may be placed on the other side of the boundaries 
to enforce the boundary conditions accordingly, as in the method of images. This results in a reduction of the overall interaction 
strength as well as the ability of each wave to counter-propagating against the mean flow. 

We believe that our observations and interpretation carries over to the non-Boussinesq Holmboe^ problem investigated 
previously^—, which is often attributed to the interaction between a Rossby wave and a gravity wave 7 ,10 ' 17 ' 24 , with wave 
speed governed by A Q and R respectively. In the piecewise-linear set up as in Holmboe’s original set up, the non-Boussinesq 
term affects the gravity waves supported on the density defect within the shear layer but not the Rossby waves, since AN 2 = 0 
at the location of the vorticity defects. Normally there is a symmetric Holmboe mode arising from the interaction between a 
leftward Rossby wave with a rightward gravity wave, together with a rightward Rossby wave and a leftward gravity wave. Non- 
Boussinesq effects will modify the gravity waves so that the interaction is no longer symmetric, and the instabilities should have 
non-zero c r , as in some previous works where the interaction was made asymmetric via other means (e.g., making the distance 
between the density and the vorticity defects asymmetric)fi r 14 i 1 7,24,30- ^34 \y e eX p ec t an analogous schematic to the one shown in 
Figure [8] should hold for the Holmboe problem. Similar effects should also be observed when smooth basic states^— 3 are con¬ 
sidered in the non-Boussinesq regime. In terms of general applicability, since large-scale stratified flows tend to be dominated 
by buoyancy effects, non-Boussinesq effects are more likely to manifest for small-scale flows. Since such instabilities has been 
observed to lead to mixin g 33 ' 36 ' 3 7 the non-Boussinesq effects on these instabilities may indirectly affect the mixing properties, 
although this possible avenue for further research is beyond the scope of this present study. 
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